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ITERATIVE  AffiTHODS  FOR  THE  SOLUTION  OF  ECONOMIC  MODELS  * 
,    by  John  Timothy  McGettigan 

^-    THE     f/ATHEI/ATICS 


Before  the  introduction  of  digital  computers,  it  was  necessary  to 
solve  economic  models  by  hand.  This  imposed  a  severe  restriction  on  the 
size  and  type  of  model  which  could  be  solved.  Although  there  are  numerous 
well  established  methods  for  solving  systems  of  linear  equations,  the 
amount  of  work  required  for  solution  increases  so  rapidly  as  the  size  of 
the  system  increases  that  hand  computation  soon  becomes  impracticable. 
Also  to  restrain  the  equations  of  the  system  to  be  linear  severely  re- 
stricts the  validity  of  an  economic  model.   Since  it  is  extremely  difficult 
to  solve  non-linear  systems  by  hand,  it  was  only  the  small,  linear,  eco- 
nomic models  which  could  be  solved  before  digital  computers  became  avail- 
able. 

The  powerful  computational  abilities  of  the  digital  computer,  in 
effect,  removed  the  restriction  of  size  from  model  building.   Using  the 
traditional  methods,  it  became  possible  to  solve  linear  systems  of  an 
extremely  high  order.   In  order  for  a  computer  to  function  it  requires  as 
algorithm,  and  these  traditional  methods  were  easily  translated  into  com- 
puter algorithms.  However,  even  with  the  use  of  a  computer,  the  solution 
of  a  non-linear  system  of  equations  is  still  quite  difficult.   This  sec- 
tion will  attempt  to  establish  the  mathematical  foundation  for  several  of 
the  algorithms  which  are  available  for  the  solution  of  such  systems.   Since 
almost  all  of  the  algorithms  available  are  extensions  of  methods  for  the 
solution  of  linear  systems,  the  mathematics  of  these  algorithms  has  been 
established  in  relation  to  their  application  to  linear  systems. 

Consider  the  general  problem  of  finding  the  solution  to  a  set  of 
linear  algebraic  equations.   Symbolically  this  involves  finding  a  vector 
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X  which  satisfies  the  equality 

AX  =  B  (1) 

where  A  is  a  non-singular  square  matrix  and  B  is  a  vector.   The  practical 
methods  for  finding  such  a  solution  can  be  divided  into  two  large  classes, 
namely,  the  exact  methods  and  the  methods  of  successive  approximations. 
In  the  first  class  are  the  simple  Gaussian  eliminations  techniques  which 
are  taught  in  high  school  algebra  classes.  Also  in  this  class  are  the 

matrix  inversion  methods  which  involve  solving  equation  (1)  for  X 

_i 
X  =  a"  B 

It  is  characteristic  of  the  exact  methods  that  they  produce  the  exact 
solution  in  a  given  number  of  operations  which  is  a  function  solely  of 
the  order  of  the  system  of  equations  and  not  of  the  matrix  of  coeffi- 
cients A. 

The  exact  methods  have  several  disadvantages  when  applied  to  rela- 
tively large  systems  of  equations.   With  these  methods  the  number  of 
operations  required  for  solution  is  roughly  proportional  to  the  cube  of 
the  number  of  equations  in  the  system:   thus  as  the  size  of  the  model  in- 
creases the  effectiveness  of  these  methods  rapidly  decreases.   Because  of 
the  limited  length  of  the  computer  registers,  there  is  a  certain  amount 
of  rounding  error  introduced  with  every  operation.  With  the  exact  methods 
these  errors  are  cumulative,  and  there  is  the  possibility  of  introducing 
significant  error.   In  general  economic  models  are  too  large  to  be  effec- 
tively handled  by  the  exact  methods:   non-linear  economic  models  which 
contained  more  than  ten  equations  could  be  more  efficiently  solved  by 
other  methods. 

Economic  models  also  have  the  characteristic  that  any  given  equation 
will  contain  relatively  few  variables.  For  a  linear  model  this  would  imply 
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that  the  matrix  of  coefficients,  A,  is  sparse  i.e.  contains  many  zeros. 
The  exact  methods  are  not  able  to  take  advantage  of  this  and  the  number 
of  operations  required  for  solution  is  not  reduced. 

The  second  class  of  methods,  the  methods  of  successive  relaxation, 
have  none  of  the  obvious  disadvantages  of  the  exact  methods,  and  therefore 
appear  to  be  better  suited  for  the  solution  of  economic  models.  The  rest 
of  this  section  will  be  devoted  to  establishing  the  mathematical  basis  for 
this  class  of  methods. 

The  general  nature  of  the  methods  of  successive  relaxation  is  con- 
tained in  the  equation 

X^-^^  =  A<">  +  f(X^^>)  (2) 

where  X  is  a  single  number,  a  function,  a  vector,  or  a  matrix,  according 

(r) 

to  context.  A    is  independent  of  X  but  not  necessarily  of  £.   The  solu- 

r+1 
tion  is  the  limit  of  the  series  of  estimates,  X   ,  as  £  goes  to  infinity. 

For  the  solution  of  a  system  of  linear  equations  the  original  vector  is 

immaterial:   if  the  process  converges,  it  will  converge  from  any  point. 

Iterative  methods  form  an  important  subset  of  the  methods  of  successive 

relaxation,  and  because  their  mathematics  is  relatively  straightforward, 

they  serve  as  a  good  introduction  to  this  area. 

Let  the  system  of  linear  equations  be  given  in  the  following  form 

Xi    =  aiiXi   +  aisxs    +  ai^x^   +  Cj 

Xn     =    SniXi     +    a„2Xs      +    a^^X^    +    Cn 

At  first  glance  such  a  representation  appears  somewhat  artificial  because 
of  the  presence  of  like  terms  on  both  sides  of  the  equality  sign.  However, 
this  arrangement  of  a  linear  system  is  the  basis  for  an  iterative  process 


-4- 
which^  although  impractical  for  applications,  is  very  convenient  for  illus- 
trating certain  properties  of  the  iterative  methods. 
Let  us  now  write  system  (1)  in  the  form 
X  «=  AX  +  C 
where  A  is  a  matrix  of  coefficients  and  C  is  a  vector  of  constant  terms. 
Starting  with  an  arbitary  initial  vector  X   ,  the  following  series  is 


formed 


X*"  ''  =  AX^  '^+  C 


X^  ^  =  AX^    ^   +  C 

If  the  sequence  X   ,  X    ,....  has  a  limit  X,  this  limit  will  be  the 

r  r+1 

solution  to  system  (2),  for  X  will  equal  X    and  this  will  produce  the 

equality  X  =  AX  +  C.  Following  through  this  recursive  function  it  is 

easily  verified  that 

X^  -*  =  A  X^  '^  +  (I  +  A  + A   )C  (3) 

Now  for  this  process  to  converge  it  is  necessary  and  sufficient  that  the 
series 

(I  +  A  + A  ~^)C  (4) 

also  converges.   It  can  be  proven  that  this  series  converges  if  and  only 
if  all  of  the  eigenvalues  of  the  matrix  A  are  of  absolute  value  less  than 
one.  And  it  is  a  necessary  and  sufficient  condition  for  the  convergence 
of  all  iterative  processes  that  the  eigenvalues  of  the  iterative  matrix 
be  of  absolute  value  less  than  one. 

For  this  method  the  solution  vector  is  the  sum  of  a  series  of  terms, 
and  referring  back  to  the  original  mathematical  description  of  the  methods 
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of  successive  approximations,  X  in  this  case  is  a  function.  The  closer 
the  largest  eigenvalue  is  to  zero  in  absolute  value,  the  faster  the  terms 
in  the  series  in  (4)  tend  to  zero  and  the  faster  is  convergence.  After  a 
few  iterations  one  should  be  able  to  estimate  the  speed  with  which  this 
sura  is  approaching  a  limit  and  hence  estimate  the  limit.   Since  the  speed 
of  convergence  is  solely  a  function  of  the  largest  eigenvalues  of  the 
iteration  matrix,  it  is  possible  to  estimate  the  solution  at  every  itera- 
tion knowing  only  the  value  of  the  largest  eigenvalue. 

Although  this  method  is  not  very  practicable  for  application  to  real 
systems,  it  shows  more  clearly  than  other  methods  that  both  the  criterion 
for  convergence  and  the  speed  of  convergence  of  the  iterative  methods  are 
a  function  of  the  largest  eigenvalue  of  the  iteration  matrix. 

The  Jacobi  iterative  methodis a  more  efficient  method,  particularly 
when  adopted  to  the  solution  of  non-linear  systems.   Let  a  system  of 
linear  algebra  equations  be  represented  in  the  matrix  notation 

AX  =  B.  (5) 

In  order  to  express  exactly  the  mathematics  of  the  Jacobi  method,  it  is 
necessary  to  decompose  the  matrix  of  coefficients.  A,  into  three  compo- 
nents, U,  L,  and  D,  which  represent  the  upper-triangular,  lower-triangular, 
and  diagonal  elements  of  A.[l]Using  this  notation  the  system  now  becomes 

(U  +  L  +  D)X  =  B. 
Taking  advantage  of  the  commutative  property  of  matrix  multiplication,  this 
equation  can  be  transformed  into 

DX  =  B  -  (U  +  L)X. 
The  Jacobi  method  approximates  this  equality  by  the  following  iteration 


-6- 

function 

r+1  r 

DX    =  B  -  (U  +  L)X 

or 

X^^  =  d'-"-  B  -  d'-'-CU  +  L)X^ 

r+1 
If  this  process  converges,  then  for  sufficiently  large  r,  X    will  not 

differ  significantly  from  X  ,  and  by  working  backwards  it  can  be  easily 

verified  that  this  vector  is  the  solution  of  the  original  system.  To  see 

what  this  process  actually  involves,  consider  the  following  representation 

of  a  linear  system 

ll    =■  aiixi    +  ai2X^..ainXn 


In    =   a^jXi    +  3^2X2  .  ..a„„Xn 

Using  this  notation  the  approximation  to  the  solution  on  the  (r+l)th 

iteration  of  the  Jacobi  method  is 

r  +  1                    2                        2 
xi        ■=     ai  p xp.        +  a^axs    -      1 


an  ai  1  an 

r  +  1  ,  r 

^n        "  a^^xi   +  apgXg      ....    -     J^ 

ann  ann  &^n 

For  the  Jacobi  scheme  the  matrix  of  iteration  is  D~  (U  +  L),  and  this 
process  will  converge  to  the  solution  vector  if  all  of  the  eigenvalues  of 
this  matrix  are  less  than  one  in  magnitude. 

The  Gauss-Seidel  and  the  Jacobi  iterative  are  closely  related. 
Decomposing  the  matrix  of  coefficients  in  the  same  manner,  equation  (5) 
is  solved  slightly  differently  to  produce  the  Gauss-Seidel  iterative  method 

■+1   _   ,,  .^-l,^    _.^  ,jr  .  ,^    ^,-1 


f^^     =      (I-L)  ^(D  +  U)  X^  +  (I  -  L)"^B 
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Using  the  more  detailed  notation,  the  solution  vector  on  the  (r+l)th 
iteration  is  approximated  by  the  vector 


Xl 

s 

ajsxs   + 

ai3  X3    ....    - 

^^ 

aix 

^11 

ai 

r+l 

S3 

r  +  l 

T 
+    %  3  ^      •  •  •  • 

-I2 

^2 

as2 

as  a 

xr' 

C3 

r  +  l 
ann 

r  +  l 
a.o  Xr)     ...... 

..    -1. 

ana 

ann 

This  process  differs  from  the  Jacobi  method  in  that  the  vector  of  approx- 
imations, X,  is  immediately  updated  with  the  new  estimates  of  the  individ- 
ual components.  Thus  for  a  system  of  order  n,  the  solution  vector  is 
updated  n  times  during  every  iteration,  whereas  for  the  Jacobi  method 
it  is  updated  only  once  for  every  iteration. 

The  matrix  of  iteration  for  the  Gauss-Seidel  method  is 

(I-  L)-^  (D  +  U) 
and  again  it  is  necessary  that  all  of  the  eigenvalues  of  this  matrix  be 

of  absolute  value  less  than  one  for  the  process  to  converge.   For  both  of 
these  methods  there  are  many  ways  in  which  the  system  of  equations  can  be 
normalized.  The  term  nonnalization  in  this  context  refers  to  which  vari- 
able a  given  equation  is  solved  for.  For  instance,  the  first  equation  may 
be  solved  for  the  tenth  variable.   If  the  matrix  of  coefficients  for  a  n 
order  linear  system  contained  no  zeros,  there  would  be  n  factorial  differ- 
ent normalizations  possible.   Since  each  normalization  produces  a  differ- 
ent matrix  of  iteration,  there  is  the  possibility  that  for  a  given  system 
the  process  will  converge  with  certain  normalizations  and  not  for  others. 
Also  since  the  solution  vector  in  the  Gauss-Seidel  method  is  continually 


updated,  the  way  in  which  the  equations  are  ordered  will  affect  the  matrix 
of  iterations.  For  the  Jacobi  method  the  ordering  is  immaterial. 

There  is  an  extensive  list  of  conditions,  which  if  satisfied,  guaran- 
tee that  the  eigenvalues  of  the  matrix  are  less  than  one  in  absolute 
value.   Unfortunately  these  are  only  sufficient  conditions  and  there  is 
a  large  class  of  matrices  for  which  the  eigenvalues  are  less  than  one  in 
absolute  value  which  do  not  satisfy  these  conditions.   Since  it  requires 
approximately  the  same  amount  of  work  to  calculate  the  largest  eigenvalue 
as  to  perform  the  entire  iteration  process,  the  determination  of  the 
largest  eigenvalue  is  an  impractical  way  of  testing  for  convergence.   Un- 
fortunately, as  yet,  there  is  no  practical  algorithm  for  normalization  or 
ordering,  and  one  can  only  proceed  on  a  hit  or  miss  basis. 

As  was  seen,  in  solving  a  system  of  linear  equations  given  in  the 
form  X  °  AX  +  F,  the  solution  vector  is  the  sum  of  a  series  of  vectors. 
Since  the  speed  with  which  this  sum  converges  to  the  solution  vector  is 
a  function  of  the  largest  eigenvalue  in  absolute  value  of  the  iteration 
matrix,  it  is  possible  at  every  step  to  estimate  the  final  sum  of  this 
series.  And  for  every  iteration  scheme,  the  speed  of  convergence  is  in- 
versely proportional  to  the  largest  eigenvalue  of  the  iteration  matrix. 
This  implies  that  when  the  largest  eigenvalue  is  close  to  one  in  absolute, 
there  is  only  a  small  step  taken  towards  the  solution  vector  during  each 
iteration.  Thus  the  speed  of  convergence  can  be  improved  by  taking  a 
larger  step  during  each  iteration  than  the  iterative  scheme  prescribes, 
and  the  closer  the  largest  eigenvalue  is  to  one  in  absolute  value,  the 
larger  the  step  should  be.  For  every  iteration  method  the  following 
approximation  holds 

1  -  e. 
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where  ej  is  the  largest  eigenvalue  in  absolute  value  of  the  iteration 
matrix.  This  approximation  can  then  be  taken  as  the  current  estimate 
of  the  solution  for  the  next  iteration.  Thus  the  number  of  steps  re- 
quired for  solution  is  significantly  reduced  when  the  value  of  the 
largest  eigenvalue  is  known. 

Besides  the  iterative  schemes,  there  is  another  class  of  methods 
which  come  under  the  heading  of  methods  of  successive  approximation. 
These  are  known  as  the  methods  of  steepest  descent,  or  the  gradient 
methods.  For  an  understanding  of  these  methods,  consider  the  relation- 
ship which  results  when  a  vector  of  initial  estimates  is  substituted 
into  a  system  of  linear  equations 

R  »  AX  -  B. 
If  X  were  the  true  solution  vector,  then  the  vector  R  would  have  only 
zeros  for  entries.  The  purpose  of  the  gradient  method  is  to  produce  a 

series  of  vectors,  the  limit  of  which  forces  R  into  the  null  vector. 
During  each  iteration  it  changes  a  given  variable  such  that  the  sum  of 
squares  of  the  components  of  the  resulting  vector  R  is  brought  to  a 
minimum.  For  instance,  if  the  residual  of  the  first  equation  is  100, 
then  by  changing  variable  one  in  this  equation  by 
100 

units,  the  residual  of  this  equation  will  become  zero.  The  general 
gradient  method  reduces  the  sum  of  squared  residuals  by  taking  the 
weighted  average  of  these  changes  for  each  variable.  The  change  in  vari- 
able i   during  iteration  is  given  by 

aiiTi  +  ajgrs +  ai„ro 

afi  +  ah +  af  „ 
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It  can  be  proven  that  for  a  system  of  linear  equations  the  gradient  method 
will  converge  from  any  starting  point;  however  the  computation  for  every 
iteration  is  lengthy  and  the  rate  of  convergence  may  often  be  slow.  The 
Jacobi  and  the  Gauss-Seidel  method  have  the  advantage  that  convergence 
can  be  very  rapid,   but     these  processes  may  not  converge. 

The  application  of  these  methods  to  non-linear  systems  is  relatively 
straightforward.   During  each  iteration  the  partial  derivatives  of  each 
variable  are  held  constant;  this  suppresses  the  non-linearities  and  the 
system  of  equations  becomes  linear.   For  example  consider  the  following 
equation 

10  =■  5x  +  wyz. 
If  the  current  estimates  for  the  solution  values  of  w,  x,  y,   and  z,  are 
respectfully  2,   A,  6,  and  8,  this  equation  would  be  forced  into  one  of 
the  following  forms 

10  =  5x  +  12z 

10  =  5x  +  16y 

10  =■  5x  +  48w 
When  the  non-linear  equations  are  linearized  in  this  manner,  it  is  possi- 
ble to  apply  the  algorithms  which  were  described  in  this  chapter.  There 
is  a  great  deal  of  mathematics  concerned  with  the  methods  of  successive 
approximations  as  applied  to  linear  systems  of  equations,  but  unfortunately 
there  has  been  little  written  about  these  methods  as  applied  to  non-linear 
systems.  To  a  great  extent  this  is  due  to  the  fact  that  the  class  of  non- 
linearities  is  so  large  that  there  is  little  which  can  be  said  about  the 
mathematics  of  such  systems.   Economic  models  contain  a  restricted  class 
of  non-linearities.  The  object  of  this  paper  was  to  apply  several  of  these 
iterative  methods,  particularly  the  Gauss-Seidel,  to  an  economic  model  to  see 
how  efficiently  they  produced  the  solution,  and  if  there  were  any  characteristics 
in  their  application  to  economic  models. 
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II.   THE  MODEL 


It  can  be  seen  by  solving  some  ordinary  systems  of  linear  equations 
that  when  the  Gauss-Seidel  method  converges,  it  converges  very  quickly. 
However,  just  as  the  mathematics  implies,  this  method  does  not  always 
produce  convergence.   On  the  other  hand,  the  gradient  method  always  con- 
verges, but  the  rate  at  which  it  approaches  the  solution  is  generally 
slow.  The  application  of  these  two  methods  to  non-linear  economic  models 
raises  some  interesting  questions.   In  what  portion  of  the  time  will  the 
Gauss-Seidel  method  converge,  and  how  fast  is  the  average  rate  of  con- 
vergence? What  is  the  best  way  to  solve  for  the  variables  and  then  order 
the  equation  in  order  to  produce  the  optimum  conditions  for  convergence 
of  the  Gauss-Seidel  method?  To  what  extent  is  convergence  conditioned 
upon  the  starting  values?  Will  the  gradient  method  always  converge  when 
applied  to  a  non- linear  system  which  has  a  real  solution? 

In  an  attempt  to  answer  these  questions,  these  methods  were  used  to 
solve  a  real  economic  model.  The  model  which  was  chosen  was  a  modifi- 
cation of  the  "Revised  Klein-Goldberger  Model. '[3"]  This  is  a  macro- 
econometric  model  which  presents  a  Keynesian  analysis  of  the  United 
States  economy.  The  equations  for  this  model  were  estimated  from  annual 
time  series  of  the  United  States  economy  for  the  period  1929  to  1962, 
with  the  years  1942  to  1945  omitted.  These  iterative  methods  were  used  on 
this  model  to  solve  a  simulation  of  the  1929  economy.   Originally  the  model 
contained  twenty-seven  equations,  but  when  used  to  simulate  the  1929  econ- 
omy, four  of  these  equations  become  recursive.   Of  the  remaining  twenty-three 
equations  twelve  are  non-linear.  With  one  exception,  all  of  the  non-linear 
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terms  are  products  which  have  the  variable  representing  "price"  as  one  of 
the  terms.  This  makes  the  system  non-linear  enough  to  give  the  methods  a 
good  test.  The  model  is  large  enough  to  place  a  premium  on  speed.  The 
difference  of  a  half  a  second  between  two  methods  when  solving  a  system 
of  five  equations  may  not  be  important,  but  as  the  size  of  the  system 
increases  the  difference  between  the  two  methods  is  likely  to  increase 
more  than  proportionally.   But  perhaps  the  most  important  characteristic 
of  the  model  is  that  it  is  a  good  representation  of  modem Keynesian 
economics;  it  is  a  representative  model. 

The  model  contains  the  following  variables: 


Consumption  of  durables,  billions  of  1954  dollars. 

Consumption  of  non-durables  and  services,  billions  of  1954 
dollars . 

Residential  construction,  billions  of  1954  dollars. 

Stock  of  inventories,  billions  of  1954  dollars. 

Imports,  billions  of  1954  dollars. 

Index  of  hours  worked  per  week,  1954  =  1.00. 

Annual  earnings,  thousands  of  dollars. 

Wages  and  salaries  and  supplements  to  wages  and  salaries, 
billions  of  1954  dollars^. 

Wage  and  salary  workers,  millions. 

Corporate  saving,  billions  of  1954  dollars. 

Corporate  profits,  billions  of  1954  dollars. 

Rental  income  and  net  interest,  billions  of  1954  dollars. 

Gross  national  product,  billions  of  1954  dollars. 

Personal  disposable  income,  billions  of  1954  dollars. 

Proprietors'  income,  billions  of  1954  dollars  plus  P  and  IVA. 
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Indirect  taxes,  billions  of  current  dollars. 

Corporate  profits  taxes,  billions  of  current  dollars. 

Personal  taxes,  billions  of  current  dollars. 

Contributions  for  social  insurance,  billions  of  current  dollars. 

Business  transfers,  billions  of  current  dollars. 

Government  transfers,  billions  of  current  dollars. 

Implicit  GNP  deflator,  1954  =  1.00. 

Inventory  valuations  adjustment,  billions  of  current  dollars. 

Investment  in  plant  and  equipment,  billions  of  1954  dollars. 

Yield  on  prime  commercial  paper,  4-6  months,  per  cent. 

Capital  consumption  allowances,  billions  of  current  dollars. 

Average  yield  on  corporate  bonds  (Moody's),  per  cent. 

Exogenous  Variables 

Net  interest  paid  by  government,  billions  of  current  dollars. 

Government  wages  and  salaries,  billions  of  1954  dollars. 

Implicit  price  deflator  for  imports,  1954  =  1.00. 

Government  employees,  millions. 

Self-employed  workers,  millions. 

Total  labor  force,  millions. 

Subsides  -  current  surplus  of  government  enterprise,  billions 
of  current  dollars. 

Dummy  variable,  0  for  1929-1946,  1  for  1947-1962. 

Average  discount  rate  at  all  Federal  Reserve  Banks,  per  cent. 

Year-end  ratio  of  member  banks'  excess  to  required  reserves. 

Government  expenditures,  billions  of  1954  dollars. 

Exports,  billions  of  1954  dollars. 

Statistical  discrepancy,  billions  of  current  dollars. 

Value  of  last  twenty  years  investment. 
Average  weekly  benefits  for  unemployed. 
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1        Cd   -  a^    Y  =   (.7  +  ag)   CcLi    -    .  7  a^    Y_i   +  83 


2        Cn  -   a^    Y  "  ag    Cn_i   + 


3        R  -  %i    Y  =  a^g   r_i   +  a^a    ILi   +  a^. 
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9  W   +    332     NW    =    W_i     +    332  (Nl     "    Ns)    +    333      (P_  ^      -    P.g)     +    33, 


10       pSc  +  IVA  -  Sge    (pPc  +  IVA   -  Tc)   =  33g    (pPc  -  Tc  -   pSc)_;^   +  a^j 
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13        X-Cd-Cn-R-H+Im=H_i+I+  G  +  E 
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14       pY  -  pX  +  pSc  +  IT  +  Tc  +  PT  +  SI   -  GT  +  IVA  =  GS   -  SD   -  D  +ID 


15        pn   -  pX  +  pW  +  pnr  +  IT  +  BT  =   -  D  +  GS   -  SD 


16        IT   -  Cg   pX  =  Ci 
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17  Tc   -  C4   pPc  =  C3 

18  (1    -  C^)    PT   -  C^    pY  =  C5 

19  SI   -  Cg    (pY  +  PT)    =  C7 

20  BT   -  ag7    pX  -  age   BT_i   +  a^^ 

21  GT  +  Cg  Nw  =  Cg   (Nl   -  Ns)  +  Cjo   WB  +  C^i 

22  pW  -  wh  Nw  =  0 

23  IVA  -  855    p  =   -  ag5    p_i   +  age 

24  I  -  a,    (X   -  Wg)_i   +  ag    r_i   +   (.95  +  a^ )    Li   +  a^^ 

25  rs  -  851  rd  +  a^^  Re.^  +853  Du  &  ag^ 

26  D  =  a^g  i:p(I  +  R)  +  a^g  Du  +  ago 

27  r  a  agg  rs  +  agg  r.^  +  837 

III.  THE  RESULTS  OF  THE  SOLUTION 


There  are  three  ways  in  which  the  Gauss-Seidel  process  can  be  varied: 
(I)  the  normalization  can  be  changed,  that  is,  different  equations  can  be 
solved  for  different  variables;  (2)  the  order  of  the  equations  can  be 
changed;  and  (3)  the  starting  estimate  can  be  varied.  From  the  mathematics 
of  this  method  as  applied  to  linear  systems,  it  is  clear  that  the  normali- 
zation is  critical;  a  system  will  converge  for  some  normalizations  and 
not  for  others.  However,  in  none  of  the  books  on  numerical  analysis  is 
there  any  mention  of  whether  the  ordering  is  a  criterion  for  convergence. 
The  model  v/as  normalized  in  approximately  ten  different  v/ays.  There  is 
only  a  limited  number  of  normalizations  possible.  Since  each  variable 
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must  be  associated  with  an  equation,  and  the  matrix  is  so  sparse  that  there 
is  not  more  than  fifteen  ways  in  which  it  can  be  renormalized.   Of  the  ten 
normalizations  tried,  only  two  were  found  to  converge.  When  the  model  was 
run  with  the  original  ordering  it  was  found  that  variable  six  was  consist- 
ently the  first  to  diverge,  eventually  pulling  the  other  variables  with 
it.  The  model  converged  when  equation  nine  was  solved  for  variable  six, 
equation  seven  was  then  solved  for  variable  nine  and  equation  six  for 
variable  seven.  With  this  normalization  approximately  eight  of  the 
twenty-seven  orderings  converged.  The  orderings  tried  were  generally 
choosen  randomly.  This  result  establishes  that  in  some  cases  the  order- 
ing is  very  critical.  The  other  ordering  which  converged  was  solving 
equation  nine  for  variable  seven  and  equation  seven  for  variable  nine. 
Thus  the  two  successful  normalizationaare  very  similar.  However,  if  only 
slight  modifications  are  made  in  these  orderings,  the  system  will  diverge. 
The  same  orderings  were  used  for  each  normalization,  and  for  this  normali- 
zation, twenty-five  of  the  twenty-seven  orderings  tried  converged. 

It  is  interesting  to  note  that  the  rate  of  convergence  for  the 
different  orderings  was  by  no  means  uniform.   Some  of  the  orderings  con- 
verged in  as  few  as  nine  iterations  while  others  required  well  over  a 
hundred.  The  criterion  for  convergence  was  rather  loose.  If  none  of  the 
solution  values  changed  by  more  than  .01  in  two  successive  iterations, 
the  process  was  said  to  have  converged.   The  initial  estimates  used  were 
the  actual  values  for  the  time  period,  but  the  process  also  converged  with 
the  zero  vector  as  the  initial  estimate. 

As  the  Gauss-Seidel  process  converges  in  each  iteration  it  will  pro- 
duce a  change  in  the  current  estimate  of  each  variable.  A  damping  factor 
was  applied  to  this  change  to  determine  whether  it  was  possible  to  force 
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convergence  upon  normalizations  and  orderings  which  normally  diverge.  With 

a  damping  factor  the  current  estimate  of  the  solution  is  given  by 

^'  =  \-l  +  ^(\  -  \-P 
where  Y'  is  the  new  estimate  of  the  vector  of  solutions,  "d"  is  the  damp- 
ing factor,  Y,  ,  is  the  estimates  solution  from  the  last  iteration,  and 
Y]^  is  the  solution  which  the  Gauss-Seidel  process  would  produce  without 
a  damping  factor.  When  the  damping  factor  is  1.  this  method  becomes 
equivalent  to  the  normal  Gauss-Seidel  process.  Damping  factors  of  .3, 
.5,  and  .7  were  tried  on  all  combinations  of  normalizations  and  order- 
ings. In  all  of  the  normalizations  for  which  the  process  had  previously 
diverged,  it  also  diverged  with  the  damping  factor.   Thus  if  the  normali- 
zation is  one  which  diverges,  applying  a  damping  factor  will  be  of  no 
benefit.   Yet  for  the  two  normalization  which  had  converged,  the  applic- 
cation  of  a  damping  factor  produced  convergence  in  every  case.  Even  for 
the  orderings  for  which  the  normal  step  diverged,  the  process  converged 
for  all  three  damping  factors.  Thus  it  appears  that  at  least  for  this 
model,  the  ordering  of  the  equations  is  not  critical  if  a  damping  factor 
is  used.   If  this  is  true  in  general,  it  in  effect  removes  one  of  the 
degrees  of  freedom  from  the  Gauss-Seidel  process.   If  the  ordering  is  not 
critical  for  convergence,  then  one  must  only  be  concerned  with  the  normali- 
zation and  the  original  estimate. 

There  is  in  general  little  lost  in  speed  when  using  a  damping  factor. 
On  the  contrary,  the  application  of  a  damping  factor  often  reduced  the 
number  of  iterations  necessary  for  convergence.  This  was  particularly 
true  when  the  normal  step  took  many  iterations  to  converge.   If  the  normal 
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step  did  not  converge  rapidly  it  was  usually  because  the  estimates  were 
oscillating  around  the  true  solution.   By  applying  a  damping  factor  in 
this  situation  the  oscillations  are  damped  out,  and  consequently  requires 
less  iterations  to  converge.  However,  when  the  normal  step  converged  in 
less  than  twenty,  the  application  of  a  damping  factor  only  prolonged  con- 
vergence.  In  all  cases  the  damping  factor  of  .7  provided  more  rapid  con- 
vergence than  the  damping  factor  of  .5,  and  similarly  the  damping  factor 
of  .5  was  more  rapid  than  .3.   In  all  cases  the  process  came  to  the  iden- 
tical solution.  The  application  of  a  damping  factor  never  forced  the 
process  to  a  false  solution. 

The  mathematics  implies  that  convergence  can  be  increased  by  taking 
larger  steps  than  the  Gauss-Seidel  method  prescribes.   Different  accelera- 
tion factors  were  applied  in  a  similar  manner  as  the  damping  factors.  A 
value  of  1.1  was  the  smallest  acceleration  factor  tried  which  produced 
convergence.   Often  when  the  normal  step  converged,  the  application  of  the 
acceleration  factor  of  1.1  made  the  process  diverge.   In  no  case  did  the 
acceleration  factor  produce  more  rapid  convergence.  Thus  as  far  as  this 
model  is  concerned,  the  acceleration  factor  never  improved  convergence. 
If  one  were  to  plot  the  number  of  iterations  necessary  for  convergence 
as  a  function  of  the  damping  or  acceleration  factor  for  each  normaliza- 
tion and  each  ordering,  the  function  in  all  cases  would  be  "U"  shaped. 
For  small  damping  factors  the  number  of  iterations  is  large.  As  the  damp- 
ing factor  increases  the  number  decreases,  and  then  sharply  increases. 
The  optimum  value  of  the  damping  function  varies  as  a  function  of  the 
normalization  and  the  ordering,  but  it  was  never  found  to  be  greater  than  1, 
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The  variables  which  were  renormalized  to  produce  convergence  were  the 
variables  for  the  index  of  hours  worked  per  week,  the  index  of  hourly  wages, 
and  the  number  of  wage  and  salary  workers.  The  two  normalizations  which 
produced  convergence  also  made  the  most  economical  sense.  The  difference 
between  the  two  successful  normalizations  were  not  great  in  terms  of 
economics.  For  the  equations  which  were  renormalized,  there  was  no  dis- 
tinct dependent  variable.   Both  of  the  successful  normalizations  made 
economic  sense.  However,  if  a  slight  modification  was  made  in  this  normali- 
zation the  process  would  diverge  immediately.  For  instance,  solving  equa- 
tion five  for  variable  thirteen  and  equation  thirteen  for  variable  five 
made  the  process  diverge.   But  this  normalization  made  little  economic 
sense,  for  it  explained  gross  national  product  as  a  function  solely  of 
imports.  Thus  for  this  model  the  only  normalizations  converged  were  the 
normalizations  which  made  economic  sense. 

From  the  results  of  solving  this  model  a  few  generalizations  can  be 
made  about  an  algorithm  for  producing  convergence  with  the  Gauss-Seidel 
iterative  method.  First  the  system  should  be  normalized  in  the  manner 
which  makes  the  most  economic  sense.  Any  ordering  may  then  be  choosen, 
but  a  damping  factor  of  approximately  .7  should  be  used.   If  the  process 
does  not  converge  then  a  new  normalization  should  be  tried.   If  a  normali- 
zation does  not  converge,  it  diverges  immediately,  so  that  it  is  immediately 
obvious  whether  the  normalization  produces  convergence.  Without  using  a 
damping  factor  one  would  have  to  attempt  many  orderings  before  it  would 
be  clear  that  it  was  the  normalization  and  not  the  ordering  which  was 
causing  the  process  to  diverge. 

If  the  model  is  to  be  solved  many  times,  it  is  profitable  to  search 
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for  an  efficient  ordering,  since  the  rate  of  convergence  between  different 

orderings  with  a  given  normalization  is  great.  The  minimum  number  of 

iterations  which  the  optimum  ordering  required  to  converge  was  the  same 

for  both  normalizations.  This  implies  that  once  a  normalization  has  been 

found  which  converges  it  is  not  necessary  to  search  for  other  normalizations 

in  order  to  speed  up  convergence.  The  rate  of  convergence  is  a  function 

chiefly  of  the  ordering.  There  is  a  pay  off  between  the  amount  of  time 

spend  searching  for  a  optimum  ordering  and  the  amount  of  time  which  this 

ordering  will  save.  The  more  times  the  model  is  to  be  solved,  the  more 

profitable  it  is  to  search  for  an  optimum  ordering.  Eventually  there  may 

be  an  algorithm  for  optimum  ordering,  but  there  does  not  seem  to  be  any 

correlation  between  the  optimum  orderings  found  in  the  solution  of  this 

model. 


-21- 


BIBLIOGRAPHY 


(1)  Faddeeva,   V.D.      Computational  Methods  in  Linear  Algebra.  New  York, 
Dover  Publications,    1959 

(2)  Khabaza,    I.M.     Numerical  Analysis.   London,    Pergamon  Press,    1965. 

(3)  Klein  Lawrence  R.      The  Keynes ian  Revolution.  MacMillan  Company,    1966. 

(4)  Ralston,  Antony.     A  First  Co-grse  in  Nujrierical  Analysis.  AfoGraw-Hill  Co. 
1965 . 

(5)  Scarborough,   J.B.     Numerical  Mathematical  Analysis.   John  Hopkins  Press, 
1962. 

(6)  Varga,   Richard,  Matrix  Iterative  Analysis.   Engelwood  Cliffs,   New  Jersey, 
Prentice-Hall,    1962. 


Date  Due 


JUL  20  *|f 

SEPgsTf 

AUG  0  4  'y-s 


AUG  2  3  198$ 


Lib-26-67 


MIT   LIBRARIES 


3    TDfiD    DD3    TSfi    771 


MIT  LIBRARIES 


TDfl 


D    D03    TE7    bSl 


MIT  LIBRARIES 


3    TDfl 


D    DD3    T56    bb4 


MIT  LIBRARIES 


3    TDfiD    003    TSfi    7m 

MIT   LIBRARIES 


3    TOfiO    003    TSfl    bE3 

MIT  LIBRARIES 


3    TDflO    003    TSfl    7b3 


MIT   LIBRARIES 


3    TOflO    003    TE7    bT3 


MIT  LIBRARIES 


1,1 
3  TOflO  D03  TSfl  fl31 


3  TOflD  003  TE7  b77 


1.1 


C^      '"^ 


vvF^-^^i 


f^^'t^-^ 


?^^ 


uZ'<&*' 


